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SECTION  I 


INTRODUCTION 


This  report  lummarizci  results  to  date  on  a  study  to  determine 
the  effectiveness  of  two  methods  which  may  be  used  to  separate  interfering 
seismic  events.  These  techniques  have  potential  value  for  countcrcvasion 
studies,  in  that  they  arc  expected  to  distinguish  underground  nuclear  explosions 
hidden  in  larger  earthquakes,  and  to  identify  sequential  underground  nuclear 
explosions  disguised  as  earthquakes. 

The  techniques  arc  applied  to  markedly  different  combinations 
of  interfering  waveforms.  The  complex  ccpstrum  technique  separates  signals 
from  the  same  azimuth,  whose  waveforms  and  thus  spectra  are  similar.  It 
does  so  by  treating  the  ertire  spectrum  of  the  data.  The  eigcnspcctrum  tech* 
nique,  by  contrast,  is  used  to  separate  signals  from  different  azimuths,  with 
greatly  different  spectra,  and  considers  the  data  over  a  narrow  frequency  band. 

The  complex  ccpstrum  technique  is  discussed  in  great  detail 
by  Schafer  (1969),  and  applied  by  him  to  removal  of  echos  in  speech.  Rooker  and 
Ong  (1972)  discuss  the  cigcnspfctrum  technique  as  applied  to  synthetic  data. 
Section  II  of  this  report  discusses  the  ccpstrum  te  chnique  as  applied  to  seismic 
data,  and  presents  results  on  artificially  constructed  interfering  ifgn.ils.  An 
analysis  of  a  possible  genuine  multiple  event  is  also  given.  Section  III  develops 
the  theory  of  the  cigcnspectrum  and  discusses  results  obtained  with  this  tech¬ 
nique  as  applied  to  interfering  chirp  waveforms  in  the  presence  of  noise. 
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SECTION  II 


COMPLEX  CEPSTRUM  TECHNIQUE 

A.  INTRODUCTION 

The  complex  cepstrum  is  a  technique  used  to  separate  a 
function  consisting  of  the  convolution  of  an  oscillatory  function  with  a  series  of 
sharp  peaks.  In  the  seismic  case,  the  oscillatory  function  is  a  waveform  ob¬ 
served  at  telescismic  distances  from  an  earthquake  or  presumed  underground 
nuclear  explosion.  The  sharply  peaked  function,  the  multipath  operator,  re¬ 
presents  either  the  multipathing  nature  of  the  earth,  or  the  fact  that  more  than 
one  seismic  event  contributed  to  the  observed  waveform. 

If  the  multipath  operator  m(t)  consists  of  a  single  delta  function 
at  time  aero,  iti  convolution  with  the  seismic  waveform  returns  just  the  ori¬ 
ginal  waveform.  If  m  consists  of  a  number  of  delta  functions  at  times  0, 

T. , . .  . ,  T  ,  with  amplitudes  1,  a . a  J-.e  convolution  yields  the  original 

In  In 

waveform  plus  the  waveform  multiplied  by  a^  and  delayed  by  T^,  plus  the  wave¬ 
form  multiplied  by  a  and  delayed  by  T_,  and  so  on. 

Cm  Cm 

Such  «  multipath  operator  might  come  about  if  the  scismi:  dis¬ 
turbance  generated  by  the  event  traveled  to  the  observation  point  by  a  number 
of  piths,  each  cf  which  attenuated  the  motion  differently.  For  surface  waves, 
continental  margins  or  mountain  chains  can  cause  multipathing.  Horizontal  re¬ 
flection  plai.es  deep  in  the  earth  can  cause  multiple  arrivals  for  bodywaves. 

Another  cause  of  multiple  arrivals  is  a  multiple  source.  Ea:th- 
quakes  are  sometimes  followed  by  aftershocks  at  quite  short  time  intervals, 
and  underground  nuclear  detonations  may  consist  of  a  number  of  events  separated 
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in  time,  for  a  variety  of  reasons.  Waves  from  such  events  need  only  follow 
one  path  to  be  observed  as  multiple  arrivals. 

The  assumption  that  the  multipath  operator  consists  of  a  series 
of  delta  functions  is  the  same  as  the  assumption  that  there  is  no  phase  distor¬ 
tion  between  the  arrivals,  and  that  they  differ  in  amplitude  but  not  in  shape. 

This  will  be  the  case  when  the  multiple  paths  followed  have  the  same  attenuation 
and  dispersion  as  a  function  of  irequency,  or  when  the  source  mechanisms 
of  the  contributing  events  are  the  same.  Obviously  these  conditions  are  never 
met  exactly,  and  the  multipath  operator  will  not  be  a  series  of  delta  functions. 
Satisfactory  results  may  still  be  achieved,  however,  if  the  later  arrivals  are 
not  too  distorted. 

D.  DESCRIP  HON  OF  THE  TECHNIQUE 

We  now  put  the  above  in  mathematical  terms.  This  has  been 
done  in  a  rigorous  fashion  by  Shafer  (1969)  and  Ulrych  (1971).  The  present 
treatment  is  intended  to  stress  the  physical  significance  of  the  various  steps. 

1.  Moduiation  of  the  Spectrum 

Suppose  that  in  the  simplest  case  the  observed  waveform  is 

given  by: 


X(t)  =  S(t)  +  a  S(t-T)  (1) 

consisting  of  a  first  arrival  S(t)  and  a  secondary  arrival  delayed  ty  T  and 
scaled  by  a.  The  second  arrival  is  a  duplicate  of  the  first  in  shape.  The 
total  signal  can  be  written  as  a  convolution: 

X(t)  =  S(t)  *  m(t)  (2) 
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where  the  multipath  operator  is  m(t)  *  g  (t)  +  a  (t-T) 
Then  the  Fourier  transform  of  the  total  signal  is: 


X{f)  =  S(f)(l  ♦  acl2"fT  ) 


Writing  the  spectrum  in  polar  form: 


X(0  *  |  S(f)  |  e1^0  /l+a2+  2a  cos2*fT  el6{0 


where  <p([)  =  tan 


-1  a  sin2"fT 


1  *  a  cos2»rfT 


(3) 


(4) 


and  *(f)  is  the  phase  spectrum  of  the  unmultipathcd  signal. 

Consider  first  the  amplitude  of  X  (f).  If  the  spectrum  of  the 
signal  alone  varies  smoothly,  the  effect  of  multiplying  it  by  >/l  +  » 2+  2a  cos  2»fT 
is  to  produce  regular  minima  at  a  frequency  spacing  of  Af  *  1/T.  For  the 
case  where  a  *  1  the  amplitude  at  the  minima  is  zero.  For  all  other  values  of 
a  the  amplitude  there  is  greater  than  zero.  We  may  say  that  the  amplitude  BpCc 
truin  is  modulated  by  the  multipath  operator. 

The  logarithm  of  the  amplitude  is  taken  before  further  pro¬ 
cessing.  The  log  amplitude  spectrum  is  then: 

In  j  S(f)  J  +  y  In  (U  a^  +  2a  cos  2*fT)  (5) 

The  modulation  is  now  additive  rather  than  multiplicative,  and  has  constant 
amplitude  across  the  frequency  band. 

Next  consider  the  phase  spectrum  for  the  special  case  where 


-1  sin  2n(T  2*fT 

<?  -  tan  : - =  — ~  \1>I 

7  1  +  cos  2fffT  2 

and  the  phase  of  the  multipath  operator  increases  linearly.  When  a  is  not  equal 
to  one  oscillations  in  the  phase  angle  appear,  ( in  addition  to  the  linear  term  of 
(6))  and,  as  with  the  amplitude  spectrum,  they  have  period  \f  c  1  /  T. 

Thus  in  both  the  amplitude  and  phase  spectrum  the  effect  of 
the  multipath  operator  is  to  introduce  a  periodic  modulation  whose  period  in 
the  frequency  domain  is  equal  to  the  inverse  of  the  delay  between  the  inter¬ 
fering  events,  and  whose  amplitude  depends  on  their  relative  amplitude. 

If  the  modulation  could  be  removed  from  the  amplitude  and 
phase  spectra,  and  those  spectra  inverse  transformed,  the  result  would  be 
the  original  waveform  S(t)  without  multipathing.  Conversely,  from  the  period 
and  amplitude  of  the  modulation  we  can  obtain  the  time  shift  and  scale  factor 
for  the  second  event.  The  ordinary  way  to  remove  unwanted  periodic  compo¬ 
nents  from  the  spectra  is  by  filtering,  and  that  is  what  is  done  in  this  ease. 

2.  Filtering  the  Spectrum 

There  arc  some  useful  analogies  between  cepstral  filtering  to 
remove  the  effect  of  the  multipath  operator,  and  filtering  a  time  scries  in  the 
frequency  domain,  if  the  roles  of  time  and  frequency  arc  interchanged.  The 
Fourier  transform  of  the  complex  logarithm  of  the  spectrum  is  called  the  com¬ 
plex  ccpstrum,  in  which  time  rather  than  fr-qucncy  is  the  independent  variable. 
Components  of  the  spectrum  periodic  in  frequency  give  rise  to  peaks  in  the 
ccpstrum  just  as  periodic  components  of  a  time  series  give  rise  to  peaks  in  the 
spectrum.  Removing  these  peaks  from  the  ccpstrum  removes  their  corres¬ 
ponding  periodic  component  from  the  spectrum. 

The  complex  ccpstrum  is  the  sum  of  the  Fourier  transform  of 
the  logarithm  of  the  amplitude  spectrum  and  the  Fourier  transform  of  the 
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phase  spectrum.  Due  to  the  symmetry  of  the  amplitude  and  phase  spectra 
the  complex  cepstrum  is  a  purely  real  function.  It  is  found  most  convenient 
to  separately  consider  these  transforms,  called  the  real  and  imaginary  cep- 
stra,  tnan  to  add  them. 

Two  general  ways  to  filter  the  cepstrum,  analogous  to  high- 
and-low  pass  filtering  and  to  bandpass  filtering,  are  used  in  this  study.  The 
low  pr.»s  filter  zeros  all  points  in  the  cepstrum  beyond  some  time  less  than 
the  cepstral  peak  it  is  trving  to  eliminate.  High  pass  filtering  retains  these 
points  while  zeroing  points  at  thort  times.  Bandp  ss  filtering  zeros  a  few 
points  centered  on  the  cepstral  peaks,  or  all  points  except  those  near  the  peaks. 
For  sho-t  period  data  it  is  found  that  the  bandpass  filter  is  more  effective  than 
the  high  or  low  p»z s  filter. 

A  simple  calculation  will  show  that  if  the  interfering  signals 
are  identical  the  cepstral  pe.k  consists  of  an  isolated  point.  Since  real  inter* 
fering  signals  are  not  precisely  identical  their  peaks  will  not  be  so  sharp  and 
it  is  desirable  to  zero  one  or  perhaps  two  points  on  each  side  of  the  peak  to 
achieve  the  best  separation. 

Peaks  in  the  cepstrum  should  be  restored  to  the  value  they  would 
have  had  if  there  were  no  multipathing.  Since  this  is  not  known,  they  are  set  to 
zero,  creating  a  disturbance  in  the  cepstrum.  The  effect  of  the  multipathing 
cannot  be  perfectly  removed  from  the  signal  for  this  reason. 

After  zeroing  the  peaks  of  the  cepstra,  the  direct  transform 
back  to  the  frequency  domain  is  taken.  Within  the  limitation  mentioned  above, 
the  spectrum  is  now  free  of  modulation,  and  represents  only  the  unmultipathed 
signal.  After  taking  the  antilogarithm  of  the  amplitude  spectrum,  an  inverse 
transform  to  the  time  domain,  called  the  shortpass  output,  recovers  the  ori¬ 
ginal  signal.  When  this  is  subtracted  from  the  input  trace,  the  remainder  is 
the  second  arrival  due  to  the  multipath  operator. 
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If  we  retain  only  the  cepstral  peaka,  and  tranaform  to  the 
requen  y  domain,  we  will  be  left  with  modulation  of  the  log  amplitude  apectrum. 
To  obtain  the  long  paaa  output,  we  firat  take  the  antilogarithm,  and  then  the 
direct  Fourier  tranaform  to  the  time  domain. 

The  mean  level  of  the  modulation  ia  zero,  becauae  the  value 
of  the  cepatrum  at  the  firat  point  waa  aet  equal  to  aero.  When  the  antilogaritnm 
ia  taken,  the  mean  value  becomea  one.  The  direct  tranaform  of  thia  conatant 
term  ia  a  point  of  amplitude  one.  There  ia  a’ao  point  in  the  long  paaa  output 
due  to  the  modulation.  The  modulation  haa  period  1/T,  ao  ita  tranaform  ia 
a  point  delayed  by  T  paat  the  point  d  ie  to  the  original  aignal.  Since  there  ia  a 
linear  term  in  the  phaae  apectrum,  due  to  the  reverae  rotation  neceaaary  to 
correctly  poaition  the  ahort  paaa  output,  both  pointa  are  ahifted  to  larger  timea 
by  a  amount  roughly  equal  to  half  the  length  of  the  original  aignal. 

3.  Interpreting  the  Reaulta 

In  practice  the  real  and  imaginary  cepatra  are  not  amooth 
functiona  aa  implied  bv  the  previoua  paragraph!.  Noiae,  preaent  to  aome  de¬ 
gree  in  all  recorda,  ia  by  ita  nature  non-repetitive,  ao  it  doea  not  produce  mo¬ 
dulation  in  the  apectra.  Rather,  it  leada  to  noiae  in  the  apectra  and  noiae  in 
the  cepatra.  Ita  net  effect  ia  of  courae  to  obacure  the  peaka  due  to  real  ar¬ 
rival  a . 

Once  the  real  and  imaginary  cepatra  have  been  calculated  and 
diaplayed,  the  analyat  muat  chooae  the  peaka  he  wiahea  to  filter  out,  which  on 
aome  ?ecorda  may  be  no  larger  than  thoae  due  to  noiae.  A  ueeful  teat  ia  to 
compare  the  real  and  imaginary  cepatra.  The  occurrence  of  peaka  at  the  aame 
time  in  both  cepatra  ia  an  indication  that  an  arrival  occura  there.  However,  it 
ia  uaually  a  matter  of  experience  and  trial  and  error  to  aucceaafully  aeparate 
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muitipathcd  signal*.  The  corrcctnca*  of  the  filler  point*  mu»t  be  judged  by 
the  quality  of  the  result*.  Some  guide*  for  such  judgement  are  given  here. 

It  ia  observed  that  when  filtering  is  done  correctly  the  short  pas* 
output  is  very  smaii  untii  the  time  cor  responding  to  the  filter  point.  If  there  .* 
no  rcai  arrival  at  this  time,  the  first  motion  wiii  not  be  at  the  filter  time  hut 
earlier. 

If  a  clean  first  motion  is  presert  one  next  iooks  for  similarity 
of  motion  between  the  output  traces.  Peaks  separated  by  the  same  time  deiay 
as  the  filter  time  should  correspond  to  one  another.  If  these  tests  arc  met,  a 

o 

negative  test  must  also  be  passed.  The  outputs  must  not  show  oscillations  180 
out  of  phase  for  any  great  part  of  the  motion.  Such  oscillations  mean  that  the 
technique  is  creating  motion  where  there  was  none  actually  recorded. 

C.  EXPERIMENTAL  RESULTS 

To  test  the  cepstrum's  ability  to  separate  identical  interfering 
kignals,  representative  waveforms  recorded  on  short  period  instruments  at 

NORSAR  were  shifted  in  time  with  respect  to  themselves,  multiplied  by  a  scale 
factor,  and  added  to  themselves.  The  multipath  operator  was  known  exactly  in 
this  ease,  and  the  performance  of  the  technique  could  be  readily  evaluated. 

Figures  II- 1  through  II -4  show  resuits  obtained  using  event 
KAZ/170/04N,  a  presumed  underground  nuclear  explosion  from  eastern  Kazakh. 
Figure  II -1  is  the  waveform,  as  recorded  with  good  signal-to-noisc  ratio  at  a 
rate  of  iO  samples  per  second.  Figure  iI-2  shows  the  result  of  shifting  the  wave¬ 
form  of  Figure  II- 1  by  .7  seconds,  multiplying  it  by  .  5,  and  adding  to  the  original 
signal.  The  result  is  a  complicated  trace  with  no  obvious  second  arrival. 

Doth  the  real  and  imaginary  ccpstra  showed  a  cicar  peak  at 
a  delay  time  of  .7  seconds  for  this  synthetic  event.  A  bandpass  filter  which 
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removed  the  points  at  .  7,  1.4...  second*  wai  applied  to  both  ccpstra,  as 
dscussed  in  Subsection  D. 

The  short  pars  output,  and  the  input  minus  the  short  pass  out¬ 
put,  are  shown  superimposed  in  F  gurc  II- 3.  The  resemblance  of  the  waveforms 
to  one  another,  and  to  the  original  tiace  of  Figure  H-l,  is  clear.  Furthermore, 
the  second  trace  is  indcco  delayed  by  .7  seconds  with  respect  to  the  fust,  and 
is  roughly  half  the  armlitud  »  of  the  first. 

The  long  pais  output,  the  estimate  of  the  multipath  operator, 

is  shown  in  Figure  II-4.  The  amplitude  of  th  ■»  s—ond  peak,  delayed  by  .7  seconds 
with  respect  to  the  first,  is  very  close  to  half  the  amplitud  of  the  first.  It  ia 
often  observed  that  the  relative  amplitudes  o.  the  long  pass  peaks  are  closer  to 
the  correct  values  than  those  in  the  short  pass  output. 

Next,  in  Figure  II- 5,  is  shown  the  waveform  of  Figure  II- 1  de¬ 
layed  by  .7  seconds,  multiplied  by  2.0,  and  added  to  itself.  Again  the  com¬ 
posite  waveform  is  complicated.  The  filtered  waveforms  are  shown  in  Figure  II- 
6,  and  the  long  past  output  in  Figure  II-7. 

The  results  in  this  case  are  not  as  clear  as  in  the  first  ex¬ 
ample,  but  the  similarity  between  the  filtered  outputs  is  still  noticeable.  The 
long  pass  output  is  again  more  clear  than  the  short  piss  output. 

Schafer  (1969;  suggests  that  the  complex  ccpstrum  technique 
is  generally  more  effective  when  the  multipath  operator  is  minimum  phase. 

When  the  second  signal  has  larger  amplitude  than  the  first,  as  in  this  example, 
the  multipath  operator  is  not  minimum  phase,  and  we  expect  less  resolution. 

Investigation  of  KAZ/ 170/04N,  and  of  KAZ/157/03N,  a  pre¬ 
sumed  underground  nuclear  explosion,  and  VAN/072/23N,  an  earthquake  near 
Vancouver  Island,  have  lead  to  the  conclusion  that  at  a  separation  of  .7  seconds, 
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FIGURE  11-7 

LONG  PASS  OUTPUT  FOR  KA7/170/04N 
SCALED  BY  2.  0  AND  SHIFTED  BY  .  7  SECONDS 
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second  events  of  relative  amplitude  from  about  2.0  down  to  a  level  where  the 
second  signal  disappears  in  the  noise  arc  separable.  The  quality  of  separation 
is  comparable  to  that  of  F’  jurcs  II- 3  and  11-6,  with  larger  second  arrivals  generally 
resulting  in  poorer  resolution  of  the  signals. 

In  the  absence  o i  an  objective  criterion  as  to  when  two  signals 
arc  resolved,  it  cannot  be  said  precisely  when  the  technique  fails.  However, 
identification  of  the  second  signal  becomes  very  dubious  when  its  amplitude  is 
more  than  twice  that  of  the  first,  at  this  reparation  in  time. 

For  these  waveforms,  even  a  tenth  of  a  second  reduction  in 
delay  time  nakes  the  signals  unrcsolvable.  We  may  conclude  that  .7  seconds 
is  the  limit  of  resolution  for  these  events. 

For  delays  of  the  order  of  1.  5  seconds  or  more  the  signals 
are  often  separable  by  eye  on  the  unprocessed  trace,  for  reasonable  scale 
factors.  For  delays  between  .  7  and  1. 5  seconds  the  upper  limit  on  the  relative 
amplitude  for  succcssfull  resolution  increases  somewhat,  but  probably  docs  not 
exceed  four.  The  lower  limit  depends  on  the  noise  as  before. 

With  these  results  in  hand  we  turn  to  a  real  event,  KAZ/115/03N. 
This  was  a  presumed  explosion,  detonated  on  April  25,  1971,  at  3  hours  32 
minutes,  58  seconds,  in  Eastern  Kazakh,  with  bodywavc  magnitude  5.9  and 
surface  wave  magnitude  4.3,  reported  by  NORSAR  (Filson  and  Dungum,  1971). 

These  magnitudes  place  it  outside  the  main  group  of  central  Asian  explosions, 
although  not  in  the  central  Asian  earthauakes. 

The  waveform  as  recorded  at  NORSAR  is  shown  in  Figure  II-8. 

Both  the  real  and  imaginary  cepstra  arc  sho.vn  in  Figure  U-9.  A  small  but  fa»rly 
ilcar  peak  is  seen  at  1.  1  seconds  time  delay  on  1<oth  traces.  Other  peaks  are 
present  on  each,  1  •  it  not  «»n  both. 
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WAVEFORM  FROM  KAZ/115/02N 
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FIGURE  n-9 

REAL  AND  IMAGINARY  CEPSTRA  FROM  KAZ/P5/03N 


A  bandpass  filter  removing  only  the  point  at  1.  1  seconds 
results  in  the  output  waveforms  shown  in  Figure  11-10.  The  lack  of  large  motion 
in  the  second  signi  1  before  1.  1  seconds  after  the  start  of  the  first  signal  is 
noticeable.  Furthermore,  the  first  four  extrema  in  the  first  waveform  occur 
just  1.  1  seconds  before  the  first  four  extrema  in  the  second  waveform,  and 
corresponding  peaks  are  in  roughly  constant  proportion  to  one  another.  At  the 
end  of  the  waveform  these  similarities  are  lost  as  noise  obscures  the  motion. 

In  the  long  pass  output,  shown  in  Figure  II-ll,  two  clear  peaks 
are  visible.  The  first  has  amplitude  one  and  the  second,  following  it  by  1.  1 
seconds,  has  amplitude  .  30.  >nother  peak  1.  1  seconds  later  has  small  am¬ 
plitude.  Thus  the  'v^dence  suggests  the  existence  of  sequential  underground 
nuclear  explosions  whose  waveforms  ar<*  separated  in  arrival  time  by  about 
l.  1  seconds,  and  whose  body  wave  amplitudes  are  about  in  the  ratio  3:  l.  Be¬ 
cause  the  motion  has  the  same  sign  for  corresponding  peaks  the  second  event  is 
not  due  to  a  surface  reflection. 

The  relative  magnitude  of  the  presumed  second  event  as  re¬ 
ported  here  is  not  sufficient  to  account  for  the  departure  of  .  3  magnitude  units 
from  the  explosion  population  reported  by  Filson  and  Bungum(l971).  However, 
the  presumed  events  need  not  have  been  detonated  at  the  same  location,  and 
differences  in  cffcciency  in  generating  body  waves  might  account  for  the  difference 
in  magnitude. 

Whether  these  waveforms  arc  due  to  two  similar  events  is  not 
a  question  the  cepstrum  technique  is  competent  to  answer.  The  claim  made  here 
is  that  the  waveforms  exhibited  in  Figure  11-10  sum  to  the  observed  v.'veform,  and 
that  the  two  individual  waveforms  are  similar.  Positive  identification  of  the 
waveforms  with  separate  presumed  explosions  is  a  step  for  the  analyst  to  take. 
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SECTION  III 

EIGENSPECTRUM  TECHNIQUE 


A.  INTRODUCTION 

About  15%  of  all  long  period  aeiarnic  events  are  "interfering", 
in  that  they  consist  of  a  mixture  of  two  different  event*  from  different  azimuth*. 
Beamsteering  van  separate  these  <v»  .its  to  some  extent,  but  if  one  is  much 
smaller  than  the  other,  it  will  be  hidden  within  the  sidelobes  of  the  main 
event.  The  eigenspectrum  technique  offers  a  possibility  of  separating  such 
events  with  greater  resolution. 

In  the  eigenspectrum  technique  the  data  matrix  is  formed  and 
its  eigenvalues  and  eigenvectors  computed.  It  is  shown  that  the  eigenvalues 
are  proportional  to  the  signal  powers  and  that  the  eigenvectors  lie  in  the 
direction  of  the  signal  vectors.  The  orthogonality  of  the  eigenvectors  makes 
possible  the  calculation  of  the  signal  powers  with  great  accuracy. 


B.  DESCRIPTION  OF  THE  TECHNIQUE 


1.  Theoretical  Basis 


For  an  array  of  sensors,  the  signal  vector  is  defined  to  be 

S  » 

where  Fjff)  is  the  Fourier  transform  at  frequency  f  of  the  motion  at  the  i  th 
site.  From  S  the  data  matrix  can  be  formed  as 


F.  (0 
•  1 


(7) 


ft  a  SS1*  * 


where  H  denotes  Hermetian  transpose. 


(8) 
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There  are  N  eigenvalues  of  this  matrix,  where  N  is  the  number  of  channels, 
and  they  are  S**S,  0  ,0  ...  Inspection  of  7  shows  that  S**S  is  N  times  the 
average  power  at  a  site.  The  eigenvector  associated  with  this  eigenvalue  is 
V. 

If  only  one  signal,  of  amplitude  A  and  wavenumber  £  is 
present  at  sites  X^,  then 

S  •  V  .  A  [jik‘*i]  ,9) 

and  V**V  is  the  power  of  this  signal.  Also,  the  eigenvector  V  is  in  the 
direction  of  the  wavevector  k.  The  more  interesting  case  is  where  two  signals 
V  and  W  are  present,  and  then 


11  r  vvH  wwH  +  vwH  *  wvH 


The  second  two  terms  are  due  to  correlation  between  the  signals.  Because  of 
their  presence,  the  eigenvalue  is 


VHV  4-  WHW  +  \VHV  +  VHW  pi) 

and  is  thus  N  times  the  total  power.  The  eigenvector  is  V+W,  and  does  not 
indicate  the  direction  of  propagation. 

The  cross  correlation  terms  may  be  eliminated  by  averaging 

the  data  matrix  U  over  frequency.  If  the  signals  are  uncor related,  the 
H  H 

average  of  VW  and  WV  will  be  zero;(this  is  the  definition  of  uncorrelated 
signals).  Under  these  conditions  the  average  data  matrix  becomes 

n.w"*W  (l2) 

The  eigenvalues  of  this  matrix  arc  not  so  easy  to  interpret 
i  >  when  the  data  matrix  represented  a  single  frequency.  We  might  hope 
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that  the  eigenvalue  of  the  average  matrix  i«  the  average  of  the  eigenvalue*  of 
the  individual  matricic*  and  thus  proportional  to  the  average  power.  However, 
this  is  so  only  when  the  eigenvectors  arc  all  parallel.  Parallel  does  not  mean 
pointing  in  the  same  direction  in  wavenumber  soacc  but: 

V(fj)  VH(f2)  x  I  VffjlM  V(f2)|  (13) 


If  we  use  (9)  for  the  signal  vector,  the  left  hand  side  of  (13)  is: 

a,VaV2>  £  x, 

i=l 

where  *  denotes  complex  conjugate 
and  the  right  hand  side  is: 

N  |  A(fj)  |  |  A(f2)| 


(H) 


05) 


and  these  quantities  are  not  equal,  although  they  approach  equality  as  f^ 
approaches 

Thus  to  obtain  a  data  matrix  which  is  free  of  cross  terms  we 
must  average  over  a  band  of  frequencies,  but  in  doing  so  the  interpretation 
of  the  eigenvalue  as  the  average  power  is  lost.  In  a  practical  situation  we 
average  over  a  band  wide  enough  to  eliminate  most  of  the  cross  correlation 
but  narrow  enough  to  not  change  the  eigenvalues  too  much  from  the  average 
power. 


The  case  where  random  noise  is  p*csent  is  easily  incorporated 
into  the  theory.  Random  noise  is  represented  as  usual  by  a  signal  vector 
whose  components  have  random  amplitudes  and  phases.  If  such  a  component 
is  n^e1^**,  then  the  j-k  th  component  of  the  data  matrix  is  n^ 

Averaging  over  frequency  reduces  those  components  for  which  j  #  k  to  zero. 
On  diagonal  components  are  n j  n.  .  The  average  over  the  frequency  hand 
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gik'cs  the  same  value  at  each  site.  Thua  the  noiae  contribution  ia  pi  where 
I  ic  the  identity  matrix.  Then  the  data  matrix  for  one  aignal  in  noiae  with  which 
it  ia  uncorrelated  ia: 

VVH  +  /?  I  (16) 

and  it  can  be  Been  that  the  eigenvalue  of  thia  matrix  la  V^V  p  ,  while  ‘he 
eigenvector  ia  V.  Thus  the  eatimate  of  the  direction  of  the  incoming  eig  ial 
ia  not  in  concept  influenced  by  noiae. 

2.  Detaila  of  the  Calculation 

The  advantage  of  the  eigenvector  technique  ia  th  it  the  eigenvaluea 
can  be  calculated  independently  of  one  another,  within  limits  of  the  calculating 
machinery.  Thua  one  aignal  ia  not  obscured  in  the  aidclobe  of  another's  beam 
pattern.  ‘ 

The  particular  eigenvalue  extraction  technique  uaed  in  thia 
study  is  called  deflation.  The  first  and  largest  eigenvector  of  a  matrix  ia 
found  by  starting  with  a  unit  trial  vector  X.  The  quantity: 

AX  *  aY  (17) 

ia  calculated,  where  Y  is  also  a  unit  vector.  Each  component  of  (X  -  Y)  ia 
compared  to  some  small  teat  value.  If  the  difference  ia  leas  than  the  teat 
value  for  all  components,  X3  Y  is  an  eigenvector  and  a  ia  its  eigenvalue.  If 
it  ia  not,  Y  is  taken  aa  the  new  trial  unit  vector,  and  the  process  is  repeated. 
Under  certain  circumatancea  the  technique  will  fail  to  find  an  eigenvalue,  but 
ordinarily  convergence  ia  achieved  in  less  than  20  iterations. 

The  data  matrix  ia  now  de  lated  by  calculating 

CtQ-  aYYH  (18) 

If  the  calculated  eigenvalue  and  eigenvector  are  close  to  their  true  values, 
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the  contribution  to  the  matrix  from  the  largest  signal  will  be  almost  completely 
removed,  and  calculation  can  proceed  in  the  lime  way  to  find  the  next  largest 
eigenvalue  and  associated  eigenvector. 

3.  Errors 

Errors  can  enter  the  calculation  at  a  number  of  points.  As 
mentioned  in  subsection  A,  the  necessity  for  averaging  the  data  matrix  over 
frequency  is  incompatible  with  the  Identification  of  the  eigenvalue  with  the 
average  power.  This  error  can  be  minimised  by  correct  choice  of  the 
frequency  averaging  band. 

Another  kind  of  error  arises  in  the  deflation  process.  The 
calculated  eigenvector  and  eigenvalue  will  not  be  exactly  equal  to  their  true 
values  due  to  limitationj  imposed  by  the  calculating  equipment.  Therefore 
the  deflation  step,  shown  in  (18),  will  leave  some  error  matrix  in  addition  to 
the  matrix  containing  the  next  largest  eigenvalue.  If  the  next  eigenvalue  is 
relatively  small,  this  error  matrix  can  significantly  distort  the  signal  power 
estimate. 

Finally,  as  a  signal  power  approaches  the  noise  level,  the 
estimate  of  the  signal  also  approaches  the  noise  power.  The  degree  to  which 
this  is  a  problem  depends  on  the  precision  with  which  we  wish  to  estimate  the 
signal  power,  as  can  be  seen  from  (16).  In  this  case,  the  error  in  the  eigen- 

spectrum  technique  is  the  same  as  in  beamsteering. 

4.  Calculations  on  Synthetic  Data 

The  work  reported  here  is  based  on  synthetic  data.  Chirp 
waveforms  were  used  to  construct  a  signal  at  each  site  oi  a  hexagonal  array. 
The  waveforms  all  had  the  form 


X(t) 


r  tin 

I2n(  t  +  2ir  Mt^  \ 

/  i  2rrt  \ 

(  1  -  cos  — r— 

\  °  T  ' 

\  T  / 

(19) 
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where  Awai  the  amplitude,  f  the  starting  frequency,  and  Af  the  frequency 
change  due  to  dispersion.  The  length  of  the  waveform  was  T,  and  the  final 
multiplicative  term  ensured  that  the  amplitude  rose  from  zero  to  A  and  then 
fell  back  to  zero. 

Each  waveform  was  Fourier  transformed,  and  the  sum  of  the 
squares  of  the  Fourier  coefficients  was  averaged  over  some  frequency  interval, 
giving  the  average  power  for  that  interval.  Then  at  each  site  the  composite 
waveform,  the  sum  of  as  many  chirp  waveforms  as  desired,  each  with  the 
proper  time  delay,  was  generated.  A  data  matrix  for  each  frequency  was 
set  up  by  calculating  the  Fourier  transform  and  forming  the  products  defined 
bv  (9).  These  data  matrices  were  averaged  to  yield  the  average  data  matrix 
defined  by  (1?).  A  random  noise  matrix  of  the  form/JI  «  as  added. 

The  largest  eigenvalue  was  extracted  from  the  data  matrix, 
and  divided  by  the  number  of  sites.  This  was  the  technique's  estimate  of  the 
power,  to  be  compared  with  the  power  found  by  Fourier  transforming  the 
original  signal.  The  eigenvector  associated  with  this  eigenvalue  was  multi¬ 
plied  by  a  look  direction  vector,  and  the  look  direction  rotated.  When  (his 
product  attained  its  largest  value  the  look  direction  coincided  with  the  signal 
direction. 

C.  EXPERIMENTAL  RESULTS 

The  signals  described  in  Subsection  D  were  generated  at  each 
site  01  a  seven  element  hexagonal  array  of  50  kilometer  apeture,  similar  to 
the  inner  ring  of  sensors  at  the  Alaska  Long  Period  Array.  The  signal  lengths 
were  different  in  order  to  ensure  that  the  signals  were  uncorrelated.  A  model 
consisting  of  a  700  second  signal  and  a  150  second  signal,  plus  random  noise, 
was  usually  assumed.  Signals  started  at  .  025  Hz  and  ended  at  .  060  Hz.  The 
azimuth  of  one  signal  relative  to  the  other  was  varied,  as  were  their  relative 
powers,  and  the  noise  power. 


111-6 


Obviously  the  technique  can  never  recover  the  original  signal 
power  and  azimuth  exactly.  It  was  decided  that  if  the  estimated  power  was 
within  .5  dR  of  the  true  power,  and  the  estimated  azimuth  within  10°  of  the 
true  azimuth,  then  the  signals  were  separated.  If  either  of  tl  esc  conditions 
were  not  met,  the  signals  were  not  resolved.  The  tolerance  on  the  power 

estimate  is  closer  than  the  azim  »th  tolerance  because  this  technique  is  ex 
pcctcd  to  be  used  to  estimate  powers  of  signals  whose  azimuths  are  already 
known  from  measurements  with  short-Dcriod  instruments. 

V/hcn  the  noise  power  was  set  equal  to  zero,  it  was  found  that 
a  signal  about  27  dR  below  a  larger  signal  could  be  resolved  over  a  wide 
angular  range.  As  the  noise  power  was  increased,  no  change  in  detectability 
was  found  until  the  noise  reached  -30  dR  relative  to  the  larger  signal.  Then 
the  estimate  of  the  smaller  signal  began  to  be  affected  by  the  noise.  As  the 
noise  power  was  raiscd»thc  level  at  which  the  smaller  signal  could  be 
detected  also  rose,  staying  3  dR  above  the  noise. 

It  was  decided  to  keep  the  noise  power  at  -30  dR  with  respect 
to  the  larger  signal.  I>owcring  the  noise  oower  would  not  increase  the  region 
of  resolution  of  the  signals,  and  raising  it  decreased  that  region. 

Figure  III- 1  shows  half  of  a  plot  of  the  region  of  detectability, 
a  small  signal  buried  in  a  large  signal,  with  the  noise  30  dB  below  the  larger 
signal.  The  other  half  of  the  plot,  from  180°  to  360°,  was  symmetrical  with  Fi¬ 
gure  III- 1 .  The  horizontal  axis  is  the  true  azimuthal  separation  between  the  two 
signals,  the  vertical  axis  is  their  relative  amplitude.  Signals  falling  within 
the  shaded  region  were  resolved  according  to  the  criteria  stated  above. 

As  the  smaller  signal  decreased  in  power  'owards  the  noise  level, 
the  eigenspertrum' s  estimate  of  power  approached  the  noise  power  asympto¬ 
tically.  As  the  azimuthal  separation  between  the  signals  decreased,  more  com¬ 
plicated  behavior  was  observed.  The  estimate  of  the  azimuth  closely  agreed 
with  the  true  azimuth  as  the  separation  decreased,  until  the  separation  reached 
about  90°.  Further  decreases  in  the  true  separation  resulted  in  errors  in  the 
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Relative  Azimuth  (degrees) 
FIGURE  IH-1 

REGION  OF  SEPARATION  FOR  INTERFERING 
CHIRP  WAVEFORMS 


estimated  seperation  of  up  to  10°.  The  estimated  separation  never  decreased  be 
low  70  ,  and  below  70°  true  difference  in  azimuth,  the  estimate  of  signal  power 
also  moved  out  of  the  tolerance  limits  for  detection.  The  cause  of  the  slight 
dtp  in  detection  level  around  180°  is  unknown. 

Behavior  at  relative  amplitudes  of  less  than  -14  dB  was  not 
investigated.  However,  beamsteering  uoulj  separate  such  signals,  and  there 
is  no  reason  to  believe  that  the  eigcnspcctrum  technique  would  not  also. 

Figure  III -2  shows  the  cigenspcc trum  for  the  largest  signal 
in  a  model  whose  parameters  arc: 


Component 

Azimuth 

Power 

Signal  1 

135° 

0  dB 

Signal  2 

295° 

-26  dB 

Noiso 

Isotropic 

-30  dB 

The  eigcnspcctrum  is  the  product  of  the  eigenvector  with  a  look  direction  vector. 
If  the  signals  are  separated  this  product  will  have  one  peak  at  the  signal  azimuth 
and  no  other  peaks.  The  eigcnspcctrum  in  Figure  IU-2  peaks  at  135°,  just  the 
signal  azimuth.  The  eigenvalue  is  almost  identical  to  the  signal  power. 

Figure  III- 3  shows  the  cigenspcct-um  for  the  second  eigenvalue. 

It  peaks  at  300°,  an  error  of  5°  in  azimuth,  but  well  within  the  criterion  for  re¬ 
solution.  The  eigenvalue  lies  25.  5  dB  below  the  first  eigenvalue,  also  within  the 
criterion  for  resolution. 

A  third  eigenvalue  was  extracted  from  the  data  matrix,  although 
only  two  signals  were  modeled.  Its  eigenvalue  lay  36.4  dB  below  the  first  igen- 
valuc,  and  its  eigcnspcctrum  is  shown  in  Figure  III -4.  The  presence  of  two  peaks 
indicates  that  the  eigenvector  is  not  of  the  form  of  a  signal  vector.  Consequently, 
this  eigenvector  and  eigenvalue  could  not  be  mistaken  for  a  signal. 
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EIGENSPECTRUM  FROM  FIRST  EIGENVECTOR 


FIGURE  IH-4 

EIGENSPECTRUM  FROM  THIRD  EIGENVECTOR 


For  comparison.  Figure  7TT-5  shows  the  beamsteer  spectrum 
for  this  signal  model.  Only  the  largest  signal  is  resolved.  The  maximum- 
likelihood  spectrum  of  Figure  111-6  resolves  both  signals,  but  is  in  error  by 
15  dB  in  its  estimate  of  the  second  signal's  power. 

From  this  comparison  it  can  be  seen  that  the  eigenspectrum 
has  substantially  greater  resolving  power  than  either  beamsteer  or  maximum- 
likelihood  spectra  for  synthetic  data. 

D.  CONCLUSIONS 

The  eigenspectrum  technique  has  the  potential  ability  to  separate 
interfering  signals  by  closely  estimating  both  their  azimuth  and  absolute  power. 
This  ability  is  baied  on  the  identification  of  the  eigenvalues  of  the  data  matrix 
with  the  signal  powers,  and  the  eigenvectors  with  the  signal  vectors.  Accurate 
algorithms  for  the  extraction  of  eigenvalues  from  matrices  can  then  be  used  to 
find  the  signal  power  and  azimuth. 

The  technique  has  been  tested  with  synthetic  data  which  simu¬ 
lates  realistic  signals  in  the  presence  of  isotropic  noise.  Over  a  wide  azi¬ 
muthal  rang'  the  eigenspectrum  technique  can  resolve  signals  down  to  27  dB 
below  another  interfering  signal.  It  is  superior  to  either  beamsteering  or  max¬ 
imum-likelihood  techniques  in  this  respect. 
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CONCLUSIONS 


Two  counter-evasion  technique*  have  been  investigated  here,  using 
synthetic  data  in  each  case.  The  complex  ccpstrum  technique  is  capable  of  se¬ 
parating  overlapping  events  which  arc  similar  to  one  another  and  whose  difference 
in  arrival  time  is  at  least  .  7  seconds,  if  the  second  signal  is  not  much  larger  than 
2.0  times  the  first.  Using  experience  gal  cd  in  processing  synthetic  combinations 
of  events,  it  is  suggested  that  the  event  KAZ/115/0*»N  may  consist  of  two  events, 
separated  in  arrival  time  by  1 .  1  seconds,  where  the  necond  event  has  about  one 
third  the  amplitude  of  tnc  first. 

The  eigcnspcctrum  technique  has  been  shown  to  have  marked 
superiority  to  both  the  beamstccr  and  maximum  likelihood  methods  in  its  ability  to 
find  small  signals  in  the  presence  of  large  ones,  for  models  based  on  chirp  wave¬ 
forms.  The  eigcnspcctrum  is  able  to  extract  a  signal  hidden  in  another  signal  of 
27  dB  greater  strength,  over  a  wide  azimuthal  range.  The  presence  of  noise  did 
not  interfere  with  the  extraction  of  one  signal  from  the  other,  unless  the  noise  was 
as  large  as  the  smaller  signal. 
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